#include "mesh.h"
#include <iostream>
#include <QDebug>


using namespace arma;

mesh::mesh(int _nElementsX, int _nElementsY, float _initialX, float _initialY , float _dx , float _dy)
{
    nElementsX = _nElementsX;
    nElementsY = _nElementsY;
    initialX   = _initialX;
    initialY   = _initialY;
    dx         = _dx;
    dy         = _dy;

    nElements = nElementsX*nElementsY;
    nNodes = (nElementsX+1)*(nElementsY+1);

    globalMatrix.clear();
    nos = new node[nNodes];

    for (int j=0 ;j < (nElementsY+1) ;j++)
    {
        for (int i=0 ; i< (nElementsX+1) ;i++)
        {
            nos[nodeIndex(i,j)].setGlobalNumeration(nodeIndex(i,j) +1);
            nos[nodeIndex(i,j)].setPositionXY(initialX+dx*i , initialY+dy*j);
//             qDebug()<< nos[nodeIndex(i,j)].getGlobalNumeration() << nos[nodeIndex(i,j)].getPositionX() << nos[nodeIndex(i,j)].getPositionY();
        }
     }

    for (int j = 0 ; j < nElementsY ; j++)
    {
        for (int i = 0 ; i < nElementsX ; i++ )
        {

            element temp;
            temp.setNumeration(elementIndex(i,j) + 1);
            temp.setNodes(&nos[nodeIndex(i,j)], &nos[nodeIndex(i+1,j)], &nos[nodeIndex(i+1,j+1)], &nos[nodeIndex(i,j+1)]);

            elementos.append(temp);
//            qDebug() << getElement(elementIndex(i,j) + 1);
        }
    }


}

node mesh:: getNode(int num) const
{
    return nos[num-1];
}

node& mesh:: getNode(int num)
{
    return nos[num-1];
}

int mesh::nodeIndex(int i,int j)
{
    return j*(nElementsX+1) + i;
}


element mesh::getElement(int num)
{
    return elementos[num-1];
}

int mesh::getNElementsX()
{
    return nElementsX;
}

int mesh::getNElementsY()
{
    return nElementsY;
}

int mesh::elementIndex(int i,int j)
{
    return j*(nElementsX) + i ;
}


void mesh::setGlobalMatrixNumerationXY()
{
    mat matAux;
    int contador=1;

    for(int i =0; i<nNodes;i++)
    {
        matAux = nos[i].getBoundaryCondition();

        if(matAux(0,0) == 0)
        {
            nos[i].setGlobalMatrixNumerationX(contador);
            contador+=1;
        }
        if(matAux(0,1) == 0)
        {
            nos[i].setGlobalMatrixNumerationY(contador);
            contador+=1;
        }     
    }
    numEq = contador-1;
}



void mesh::setGlobalMatrix()
{
    int i=0,j=0;
    mat localMatrix;
    globalMatrix = zeros(numEq, numEq);
    for(int e=1; e < nElements + 1 ; e++)
    {
        elementos[e-1].setLocalMatrix();
        localMatrix = getElement(e).getLocalMatrix();

        for(int a=1; a<5 ;a++)
        {
            for(int b=1; b<5 ; b++)
            {
                i = getElement(e).getNode(a)->getGlobalMatrixNumerationX();
                j = getElement(e).getNode(b)->getGlobalMatrixNumerationX();
                if(i!=0 && j!=0)
                {
                    globalMatrix(i-1,j-1) += localMatrix(2*(a-1),2*(b-1));
                }

                i = getElement(e).getNode(a)->getGlobalMatrixNumerationY();
                j = getElement(e).getNode(b)->getGlobalMatrixNumerationY();
                if(i!=0 && j!=0)
                {
                    globalMatrix(i-1,j-1) += localMatrix(2*(a-1)+1,2*(b-1)+1);
                }

                i = getElement(e).getNode(a)->getGlobalMatrixNumerationX();
                j = getElement(e).getNode(b)->getGlobalMatrixNumerationY();
                if(i!=0 && j!=0)
                {
                    globalMatrix(i-1,j-1) += localMatrix(2*(a-1)+1,2*(b-1));
                }

                i = getElement(e).getNode(a)->getGlobalMatrixNumerationY();
                j = getElement(e).getNode(b)->getGlobalMatrixNumerationX();
                if(i!=0 && j!=0)
                {
                    globalMatrix(i-1,j-1) += localMatrix(2*(a-1)+1,2*(b-1));
                }
            }
        }
    }
}



void mesh::setGlobalForce()
{
    mat localForce;

    globalForce = zeros(numEq,1);

    for(int e=0;e<nElements ;e++)
    {
        elementos[e].setLocalForce();
        localForce = getElement(e+1).getLocalForce();

        for(int a = 0; a<4;a++)
        {
            int i = getElement(e+1).getNode(a+1)->getGlobalMatrixNumerationX();
            int j = getElement(e+1).getNode(a+1)->getGlobalMatrixNumerationY();
            if(i!=0 && j!=0)
            {
                globalForce(i-1,0)+=localForce(2*a ,0);
                globalForce(j-1,0)+=localForce(2*a + 1,0);
            }
         }
    }

}

mat mesh::getGlobalMatrix()
{
    return globalMatrix;
}

mat mesh::getGlobalForce()
{
    return globalForce;
}


void mesh::setNodeDeslocation(vec deslocation)
{

    int contador=0;
    for(int i=0; i <nNodes ; i++)
    {
        mat aux=nos[i].getBoundaryCondition();
        if(aux(0,0)!=0 || aux(0,1)!=0)
        {
            nos[i].setHorizontalDeslocation(aux(1,0));
            nos[i].setVerticalDeslocation(aux(1,1));

        }else
        {
            nos[i].setHorizontalDeslocation(deslocation(contador));
            nos[i].setVerticalDeslocation(deslocation(++contador));
            contador++;
        }
    }
}

void mesh::setNodeDeformation()
{
    mat dPhi;
    mat aux;
    float e11,e12,e21,e22;

    for(int i=0;i<nElements;i++)
    {

        dPhi = elementos[i].getDphi();
        for(int j=0;j<4;j++)
        {
            aux=zeros(2,2);
            e11=0.5*(elementos[i].getNode(j+1)->getHorizontalDeslocation()*dPhi(j,0) + elementos[i].getNode(j+1)->getHorizontalDeslocation()*dPhi(j,0));
            e12=0.5*(elementos[i].getNode(j+1)->getHorizontalDeslocation()*dPhi(j,1) + elementos[i].getNode(j+1)->getVerticalDeslocation()*dPhi(j,0));
            e21=0.5*(elementos[i].getNode(j+1)->getVerticalDeslocation()*dPhi(j,0) + elementos[i].getNode(j+1)->getHorizontalDeslocation()*dPhi(j,1));
            e22=0.5*(elementos[i].getNode(j+1)->getVerticalDeslocation()*dPhi(j,1) + elementos[i].getNode(j+1)->getVerticalDeslocation()*dPhi(j,1));


            aux(0,0) = e11;
            aux(0,1) = e12;
            aux(1,0) = e21;
            aux(0,1) = e22;

            elementos[i].getNode(j+1)->setDeformation(aux);
        }

    }
}

void mesh::setNodeDeformationAnal()
{
    mat aux = zeros(2,2);
    float du1dx1, du1dx2, du2dx1,du2dx2;
    for(int i=0;i<nNodes;i++)
    {
        du1dx1 = -(sin(datum::pi*nos[i].getPositionX())*cos(datum::pi*nos[i].getPositionY()))/datum::pi;
        du1dx2 = -(cos(datum::pi*nos[i].getPositionX())*sin(datum::pi*nos[i].getPositionY()))/datum::pi;
        du2dx1 =  (cos(datum::pi*nos[i].getPositionX())*sin(datum::pi*nos[i].getPositionY()))/datum::pi;
        du2dx2 =  (sin(datum::pi*nos[i].getPositionX())*cos(datum::pi*nos[i].getPositionY()))/datum::pi;




        aux(0,0) = 0.5*(du1dx1 + du1dx1);
        aux(0,1) = 0.5*(du1dx2 + du2dx1);
        aux(1,0) = 0.5*(du2dx1 + du1dx2);
        aux(0,1) = 0.5*(du2dx2 + du2dx2);

        nos[i].setDeformationAnal(aux);
    }
}
